function  [H, h] = LPWIN(Wc,N,WT,P,specs)
%
%  [H, h] = LPWIN(Wc,N,WT,P,specs)
%
%  Design of Lowpass FIR Filters using WINDOWS.  This function does all
%  the number crunching involved in one iteration of filter design.
%
%  Wc = cutoff frequency of ideal lowpass to be truncated
%  N  = length of filter desired
%  WT = window type to use:
%       'rect'  'tria'  'hann'  'hamm'  'blac'
%  P  = optional variable, if nonzero a plot of the filter is displayed
%  specs = optional vector containing the specs to be met by the filter
%        = [ K1 w1 K2 w2 ]
%                            w1         w2          pi
%        |H|  0 dB |---------|----------|-----------|
%            K1 dB |---------           |
%                  |         |          |
%                  |         |          |
%            K2 dB |         |           -----------
%  H  = FFT of the impulse response of the filter
%  h  = impulse response of the filter

if ( Wc<=0 | Wc>=pi | N < 3 | WT > 5)
        error('Invalid input parameters specified')
end
if nargin < 4, P = 0; end
if nargin == 5, S = 1;end
if S&length(specs)~=4,
        error('variable specs incorectly given')
end
nm1 = N - 1;
nm1d2 = nm1/2;
n = 0: nm1;

% calculate h(n) = hd(n) * w(n)
h =  (sin(Wc*(n - nm1d2)) ./ (pi * (n - nm1d2))) .* window(N, WT);
if rem(N,2), h((N+1)/2)= Wc/pi;end      % get rid of NaN in data for N odd

m = max([ 128 (2^(ceil(log(N)/log(2)))) ]);     % find length of FFT

H = fft([h zeros(1,(m-N))]);           % get fast fourier transform of h(n)

if P,
        fftplot(H);  % plot 20*log(H(ejw))

        if S,

             hold
             a=axis;
             axis(a);
             plot([0 specs(2)/pi specs(2)/pi],[specs(1) specs(1) a(3)],':')
             plot([specs(4)/pi specs(4)/pi a(2)],...
[a(4) specs(3) specs(3)],':')
             pause
             axis
             hold

        end
end
